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abstract 


J.  M.  Richardson,  A.  B.  Arons,  and  R.  R.  Halverson  (Journal  of 
Chemical  Physics.  Vol.  15,  1947)  developed  a  calculation  procedure  for 
determining  the  hydrodynamic  properties  of  sea  water  at  the  front  of  a 
shock  wave.  The  procedure  has  been  programmed  for  the  HP41CV,  which  is  a 
hand-held  programmable  calculator .  The  program,  which  uses  374  lines 
of  code,  reproduces  the  values  for  a  shock  wave  as  tabulated  by  Richardson, 
et  al.  The  advantage  of  the  HP41CV  program  is  that  properties  can  be 
calculated  without  use  of  tables.  Copies  of  the  four  magnetic  cards  which 
have  the  program  stored  can  be  obtained  on  request. 
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I.  INTRODUCTION 


A  basic  building  block  for  describing  underwater  explosions  is  the 
hydrodynamic  properties  across  a  shock  wave  in  sea  water.  Underwater 
explosions  were  very  extensively  studied  in  World  War  II.  Results  from 
the  World  War  II  studies  are  reported  in  Underwater  Explosion  Research  [1]. 
Richardson,  Arons,  and  Halverson  [2]  discuss  the  calculation  of  hydro- 
dynamic  properties  of  saa  water  at  the  front  of  a  shock  wave.  Background 
information  on  the  calculations  appears  in  Chapter  2  of  Cole  [3]. 

Richardson,  et  al. ,  [2]  used  graphical  techniques  to  Interpolate 
thermodynamic  data;  the  methods  were  crude  and  extremely  tedious.  With  the 
advent  of  the  modern  programmable  hand-held  calculator,  the  results 
obtained  by  Richardson,  et  al..  [2]  can  be  obtained  with  ease.  This  report 
discusses  a  program  for  the  HP41CV,  which  uses  the  data  of  Richardson, 
et  al..  [2]  to  calculate  shock  wave  velocity,  water  velocity,  specific 
volume,  etc.  Reference  [2]  is  Appendix  A. 

Since  World  War  II,  research  has  continued  in  the  area  of  underwater 
explosions.  A  more  recent  survey  by  Professor  Holt  [4]  appears  in  the 
Annual  Reviews  of  Fluid  Mechanics.  Holt  [4]  discusses  underwater  nuclear 
explosions  as  well  as  chemical  explosions. 

In  passing,  an  interesting  observation  can  be  made  concerning  the  level 
of  talent  working  on  underwater  explosions  ip  World  War  II.  The  names 
include  G.  I.  Taylor,  P.  Bridgman,  H.  Beths,  J.  von  Neumann,  and  L.  I.  Sedov. 


II.  OUTLINE  OF  CALCULATION  PROCEDURE 


A.  Iteration  of  Properties  with  Rankine  Hugonioc  Relations 

The  computer  program  has  an  Iteration  which  is  now  discussed.  The 
enthalpy  jump  across  the  shock  wave  is  given  by 


AH  ■  -|(p  -  pft)  (v  +  vft) 


(2.3) 


’Q/ 

The  same  equation  numbers  and  the  same  symbols  used  by  Richardson,  et  al. ,  [2] 


are  used  here.  Since  p  »  pQ,  pQ  can  be  neglected.  Also  AH  is  given  by 


AH  -  ui  +  h  (2.20) 

where  u>  is  the  undissipated  enthalpy  and  h  is  the  dissipated  enthalpy.  The 
dissipated  enthalpy  can  be  calculated  two  ways.  First, 

hg  •  AH  -  w 

Second, 


A  H  «  h_  -  f  L  C 
P  T  T„ 


0  p 


dT 


(2.9) 


The  symbol  hg,  which  does  not  appear  in  reference  [2],  denotes  the  value 
calculated  from  AH  -  The  symbol  h^  is  the  value  calculated  from 
equation  (2.9).  Temperatures  Tq  and  T^  are  defined  in  reference  [2].  The 
undissipatad  enthalpy  is  a  function  of  specific  volume,  v,  behind  the  shock 
wave;  the  function  is 


n-1 


C  V 

“  *7^1  -  1]  (2.19) 

The  Talt  equation  of  state  relating  specific  volume  and  pressure  for  water  is 


P  ■  A[S]!(v,/v)-  1] 


where  A[S]  is  a  function  only  of  entropy.  Further 


A[S]  ■  c^/nv^ 


(2.14) 
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Combining  equations  (2.3),  (2.19),  and  (2.14)  yields 
2  2 

H  m  [(vi/v>n  "  13  (v  +  V  ■  T^I[(vi/v)a"1  “  11 


(2.21) 


Equation  (2.21)  indicates  that  is  a  function  of  v  and  v^.  The  specific 
volume,  v^,  is  a  function  of  T^.  By  varying  and  v,  the  value  of  hg 
varies. 

When  conditions  across  Che  shock  wave  are  correct, 

tijCl^)  -  h^v.Tj.) 


One  varies  v  and  T^  until  equality  is  achieved;  since  v  •  v(p,T1),  only 
needs  to  be  varied. 

This  is  a  very  brief  description  of  the  iteration.  The  reader  should 
refer  to  reference  [2]  for  more  informs cion. 

B.  Temperature  Behind  Shock  Wave 

The  temperature  increase  behind  the  shock  wave  was  calculated  using 
equation  (II~4)  of  reference  [2].  Although  the  equation  is  somewhat 
lengthy,  the  calculation  is  straightforward  except  for  two  quantities: 
B'CCq)  and  3q.  A  cubic  polynomial  is  given  for  B(tg);  this  polynomial  was 
differentiated  to  obtain 


The  quantity  3q  is  the  coefficient  of  thermal  expansion  times  specific 
volume  and  was  evaluated  using 


vhers  v  is  spsclflc  volume.  The  symbol  t  is  used  for  °C  and  T  for  °K. 
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III.  DESCRIPTION  OF  COMPUTER  PROGRAM 
A.  Computer  Program  Flov  Chare 

A  flov  chare  of  eha  program  la  given  in  Figure  1.  An  ieeraelon  loop 

occurs  for  asking  •  HH  • 

for  acceptable  difference  baeveen  HH  and  HT.  Since  HH  and  HT  have  units 

of  Joule/kg,  "LOOP"  also  has  units.  The  iteration  variable  is  t^,  which 

is  the  adiabatic  temperature  in  °C. 

When  an  acceptable  value  for  e^  has  been  found,  the  program  auto- 

3 

matically  goes  to  subroutine  SHOCK.  SHOCK  calculates  v,  cm  /gm;  u,  m/sec; 

U,  m/sec;  and  c,  m/sec. 

Following  SHOCK,  the  program  moves  to  subroutine  DEL  T,  which  calculates 
the  difference  in  eemparaeure  from  behind  to  front  of  shock  wave.  A  variety 
of  subroutines  are  used  by  DEL  T  including  CP,  3 TP RIME,  Gl,  and  D2. 

3.  Assignment  of  Storage  Registers 

For  the  program,  34  storage  registers  are  used.  The  assignment 
of  variables  to  the  registers  is  given  in  Table  I. 

C.  Program  Listing 

Appendix  B  is  a  listing  of  the  program. 

17.  OPERATION  OF  PROGRAM 

The  program  has  been  recorded  on  magnetic  cards  and  can  be  inserted 
into  the  HP41CV  using  a  HP82104A  card  reader.  The  inputs  to  the  program 
will  now  be  described. 

1.  Turn  on  the  HP41CV  and  activate  USER. 

2.  Assign  tht  program  to  a  key  on  the  key  board  by  praesing 

Q  ALPHA  HH  ALPHA,  SIN 


hj  »  HT.  Tha  input  T.OOP  is  the  criterion 


4 


!  INPUT  LOOP  1 


INPUT  Nil 


CALCULATE  j 
a+1  .  1  11 
n-1*  1  "  a’  ~  n 


INPUT  INITIAL 
TEMPERATURE 


IN°C  ; 


INPUT  FINAL 
PRESSURE  IN  KILOBARS 


XEQ  Y 

CALCULATE  (v^v)0 
FOR  USE  IN  EO.  (2.21) 


XEQ  VW 
CALCULATE  v. 


CALCULATE  hg 


I 


DISPLAY  HH,  ST 
AND  HH-HT  » 


ABS (HH-HT)  <  LOOP? 
+ 

YES 


YES 


GO  TO  "SHOCK"  1 

CALCULATE  SPECIFIC  J 
VOLUME  v  BEHIND  j 
SHOCK  wave  : 


CALCULATE  SPEED 
OF  WATER  BEHIND 

SHOCK  WAVE,  u 

- n5 — - 


’I 


CALCULATE  SPEED 


- I  " — 

OF  SHOCK  WAVE.  U 

CALCULATE  vQ 
XEQ  VW 

_ _ 

. 

CALCULATE  SPEED 
OF  SOUND  3EHIND 
SHOCK  WAVE,  c 

I— ! - ‘ - 1 

XEQ  DEL  T 


XEQ  CP 
CALCULATE  C 


jlJ 


XEQ  BT  PRIME 
CALCULATE  dB/dt 


. T" ...  I 

XEQ  G1 

XEQ  SPEED 

calculate  g 

CALCULATE  c*/n 

_ i _ 

XEQ  D2 

1 

CALCULATE  D 

XEQ  HT 

.  i  . 

CALCULATE 

!  CALCULATE  AT  1 

"1  ~ 

DISPLAY  AT  ! 
I  END  1 


Figure  1.  Computer  Flow  Chert. 
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Table  I.  Storage  Regis cers 


Register  Symbol 

Definition 

Units 

Programs 

01 

C1 

adiabatic  temperature 

°C 

HT,  VW,  HH, 
SPEED 

02 

(t-55)2 

quantity  in  polynomial 

°c 

BT 

03 

(t-25) 

quantity  in  polynomial 

°c 

VW 

04 

(t-55) 

quantity  in  polynomial 

°c 

BT 

05 

C0 

temperature  in  front  of 
shock  wave 

°c 

HH 

06 

“t 

HT,  dissipated  enthalpy 

Joule/kg 

HT 

07 

HE,  dissipated  enthalpy 

Joule/kg 

HH 

08 

B(t) 

quantity  in  Tait  equation 
of  state 

kilobars 

BT 

09 

8(t) 

quantity  in  Tait  equation 
of  state 

N/m2 

BT 

10 

t-55 

quantity  in  polynomial 

°C 

BT  PRIME 

11 

A[S] 

quantity  in  Talc  equation 
of  state 

N/m2 

SPEED 

12 

P 

pressure  behind  shock  wave 

kilobars 

Y 

13 

V 

specific  volume 

3 

cm  /gm 

SHOCK,  VW 

14 

u 

velocity  of  water  behind 
shock  wave 

m/sec 

SHOCK 

15 

p 

pressure  behind  shock  wave 

N/m2 

SHOCK 

16 

0 

velocity  of  shock  wave 

m/sec 

SHOCK 

17 

ci 

speed  of  sound  at  adiabatic 
temperature 

m/sec 

SHOCK 

18 

e 

speed  of  sound  behind  shock 

wave 

m/sec 

SHOCK 

19 

S0 

coefficient  of  thermal 
expansion  times  specific 
volume  Vq 

m3/kg°K 

BETA 

Table  1  Continued.  Storage  Registers 


Register  Symbol 


Definition  Units  Programs 


20 

o 

H  N) 

P 

quantity  in  equation 
(2.21) 

mVsec2 

SPEED 

21 

(04-1)  /  (n-1) 

function  of  n 

- 

HH 

22 

1-1/n 

function  of  n 

- 

HH 

23 

-1/n 

function  of  n 

- 

HH 

24 

dB/dt 

derivative  of  3(t) 

S/a2  °K 

BT  PRIME 

25 

y 

(v^/v)n  ratio  of 
specific  volumes 

- 

HH,  Y 

26 

V1 

specific  volume  at 
adiabatic  temperature 

ea3/gm 

HH,  VW 

27 

vo 

specific  volume  in 
front  of  shock  wave 

3. 

cm  /gm 

HH 

28 

n 

exponent  in  Tait  equa¬ 
tion  of  state 

- 

HH 

29 

C 

P 

heat  capacity  at 
constant  pressure 

goule/kg 

X 

CP 

30 

G 

quantity  in  equation 
(II-4) 

°K 

G2,  DEL  T 

31 

D 

quantity  in  equation 
(II-4) 

- 

D2,  DEL  T 

32 

1+p/B 

quantity  in  Talt  equa¬ 
tion  of  state 

- 

HH 

33 

(1+p/B) 1/n 

quantity  in  Tait  equa¬ 
tion  of  state 

- 

HH 

34 

LOOP 

criterion  for  acceptable 
difference  between  HE 
and  HI;  typical  value 

1.0 


Joule/kg  HH 


r 


la  this  example,  the  program  has  been  assigned  co  the  SIN  key 
with  blue  H. 

3.  Press  SIN  and  observe  LOOP?.  LOOP  is  the  value  of  HH-HT  which 
is  acceptable.  When  ABS (HH-HT)  is  less  than  LOOP,  the  program 
is  out  of  the  iteration  loop.  A  typical  value  for  LOOP 

could  be  1.0  Joule/kg.  Values  for  a  sample  case  are  given  in 
parentheses.  (1.0) 

4.  Press  R/S  and  observe  N?.  N  is  the  exponent  in  Tait  equation 
of  state.  Holt  (4]  uses  7.0.  Richardson,  et  al..  [2]  use  7.15. 
Insert  a  value  for  N.  (7.15) 

5.  Press  R/S  and  observe  TEMP  0?.  This  is  the  temperature  in  front 
of  the  shock  wave  in  °C.  Insert  a  value.  (0) 

6.  Press  R/S  and  observe  PRESSURE?.  This  is  the  pressure  behind  the 
shock  wave  in  kilobars.  Insert  a  value.  (80) 

7.  Press  R/S  and  observe  TEMP  1?.  This  is  the  adiabatic  temperature 

discussed  in  reference  [2].  Insert  a  value  in  °C.  (180) 

3.  Press  R/S  and  observe  crows  foot  moving  across  the  display.  In 
10  seconds,  a  value  for  HH  appears.  (520,245) 

9.  Press  R/S  and  observe  value  for  HT.  (724,255) 

10.  Ptsss  R/S  and  observe  value  for  HH-HT.  (-204,010) 

Since  Che  value  for  HH-HT  exceeds  LOOP,  a  new  value  for 
t^  •  TEMP  1  is  needed. 

11.  Press  R/S  and  observe  TEMP  1?.  To  decide  what  value  of  TEMP  1 
to  insert,  look  at  equation  (2.9);  h^  increases  as  T^  or  TEMP  1 
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increases.  Although  not  readily  apparent  from  equation  (2.21), 
hg  decreases  as  TEMP  1  increases.  For  the  example,  a  small 
table  can  be  made  as  follows: 


TEMP  1 

HH 

HT 

HH-HT 

iao 

520,245 

724,255 

-  204,010 

160 

569,362 

642,846 

-  72,984 

149.2 

591,574 

598,988 

-  7,414 

148.0 

593,321 

594,120 

-  299 

147.75 

594,285 

593,106 

1,130 

147.949 

593,916 

593,913 

3 

147.9496 

593,915 

593,915 

-  0 

12.  Press  R/S  and  observe  the  following: 

▼  *  0.6576  (ca'Vgm) 

WATER  U  -  1,638  (a/sec) 

SHOCK  U  •  4,891  (a/sec) 

WAVE  C  -  6,264  (a/sec) 

13.  Press  R/S  and  observe  the  following: 

DEL  T  -  371.24°K 

Hence  the  temperature  behind  the  shock  wave  is 
T  *  TQ  +  AT  »  273.16  +  371.24  -  644. 4°K 
The  calculation  of  shock  properties  has  now  been  completed.  To  calculate 
values  for  another  set  of  tg  and  p,  press  SIN  key  with  blue  H.  The 
computer  is  in  USER  mode. 
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V.  SAMPLE  RESULTS 


A  comparison  can  be  made  between  results  calculated  with  the  program 
and  values  tabulated  by  Richardson,  et  al. ,  [2].  Table  II  is  a  comparison. 
Most  values  are  within  1/2  per  cent. 

Since  Richardson,  et  al..  £2]  used  graphical  techniques  for  interpola¬ 
tion,  the  calculator  values  probably  are  more  accurate.  However,  the 
thermodynamic  data  for  water  have  been  represented  by  polynomials.  Errors 
in  curve  fitting  are  undoubtedly  a  few  per  cent.  Hence,  one  would  expect 
the  values  predicted  here  to  be  correct  within  a  few  per  cent. 

VI.  MAGNETIC  CARDS 

Copies  of  the  magnetic  cards  containing  this  program  are  available 

on  request  from  the  following  address: 

Distinguished  Professor  Allen  E.  Fuhs 
Department  of  Aeronautics,  Code  67 
Naval  Postgraduate  School 
Monterey,  CA  93940 

(408J-646-2948 

AV-878-2948 

VII.  CALCULATION  OP  WATER  PROPEP.TIES 
DURING  ISENTROPIC  EXPANSION  FROM  SHOCKED  CONDITIONS 

In  order  to  calculate  the  complete  flow  field  due  to  a  shaped  charge 
jet  penetrating  water,  the  calculation  of  temperature,  internal  energy,  and 
enthalpy  during  iseatropic  expansion  is  necessary.  Appendix  C  discusses  two 
subroutines  which  were  developed  after  the  main  program  was  written.  The 
subroutines  calculate  internal  energy  and  enthalpy,  subroutine  E-&AR,  and 
temperature,  subroutine  I SEN  T. 
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Table  II.  Properties  of  Sea  Water  at  a  Shock  Front* 


APPENDIX  A 


HYDRODYNAMIC  PROPERTIES  OF  SEA  WATER  AT  THE  FRONT  OF  A  SHOCK  WAVE* 


*Thia  1*  reference  [2]  reproduced  from  reference  [1].  The  technical  paper 
la  provided  aa  a  convenience  to  the  reader. 
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Hydrodynamic  Properties  of  Sea  Water  at  the  Front  of  a  Shock  Wave* 

J.  M.  Richauson** 

Saktr  Laboratory,  Cor noil  i  narrmy,  Ithaca.  S.  Y. 

AND 

A.  B.  Axons***  and  R.  R.  Halvbxson**** 

Undtrmaur  Expiosms  Rttoarck  Laboratory,  Woods  Holt  Ocoanopapkic  Institution,  Woods  Holt,  MassackuttUs 

(Received  June  2,  1947) 

The  Ranlcine-Hugonioc  relations  have  been  applied  to  appropriate  equation -of-stace  data 
in  order  to  calculate  the  propagation  velocity,  particle  velocity,  enthalpy  increment.  Riemann 
(unction,  etc.  at  shock  fronts  of  various  amplitudes  in  sea  water.  One  set  of  tables  provides 
values  over  a  wide  pressure  range  (up  to  about  30  Icilobarsi  and  is  principally  intended  for  use 
in  conjunction  with  theories  of  propagation  of  shock  waves  originated  by  underwater  ex¬ 
plosions.  A  second  set  of  tables  contains  values  which  are  closely  spaced  up  to  pressures  of  14 
Idiobara.  These  are  calculated  with  somewhat  greater  precision  and  are  intended  for  use  in  con¬ 
nection  with  experimental  measurements  of  particle  and  propagation  velocities,  etc. 


L  IHTRODOCTIOIT 

T  haa  long  been  recognized  that  the  velocity 
of  propagation  of  sound  waves  of  finite 
amplitude  in  a  fluid  medium  is  a  function  of  the 
pressure  in  the  wave.  Lamb1  ascribes  the  early 


*  The  work  described  in  this  report  was  performed 
under  .Vationnl  Defense  Research  Committee  Contracts 
OEMsr-121  with  Cornell  University  and  OEMir-569  with 
the  Woods  Hois  Oceanographic  institution. 

“  Present  address:  Mil  Telephone  Laboratories.  Inc., 
Murray  Hill.  .V.  J. 

*—  Present  address:  Department  of  Phvsici,  Stevens 
Institute  of  Technology,  Hoboken.  N.  J. 

,,M  Dectftatd* 

1  H.  Lamb.  Hydrodynamics  Cambridge  University  Press, 
London,  1932)  6th  Ed.,  p.  481. 


development  of  the  theory  to  independent  inves¬ 
tigations  of  Eamahaw  and  Riemann.  Qualita¬ 
tively  this  work  indicated  that,  since  the  higher 
pressure  portions  of  a  wave  travel  with  greater 
velocity,  an  arbitrarily-shaped  pressure  pulse  of 
finite  amplitude  must,  during  propagation,  alter 
its  shape  in  such  a  manner  as  to  build  up  into  a 
shock  front.  By  applying  the  taws  of  conserva¬ 
tion  of  mass,  energy,  and  momentum  to  the 
transfer  of  matter  across  the  shock  front. 
Rankine  and  Hugoniot  obtained  a  set  of  three 
relations  among  the  five  variables:  pressure, 
density,  particle  velocity  (u),  shock  front 
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Fig.  1.  Adiabatic  and  Hugoniot  contour*  in  ?—  T  plane. 

velocity  (If),  and  enthalpy  increment  <lH). 
These  relations,  when  applied  to  data  on  the 
equation-of-state  and  specific  heat,  make  it  pos¬ 
sible  to  calculate  u,  U ,  and  \H  and  to  evaluate 
certain  other  functions  applicable  to  the  theory 
of  the  formation  and  propagation  of  shock  waves 
originated  by  explosions.1 

Precise  knowledge  of  u  and  V  also  makes  it 
possible  to  calculate  shock  wave  pressures  in 
cases  where  the  particle  velocity  or  propagation 
velocity  can  be  measured.  The  purpose  of  the 
calculations  described  below  was  to  apply  the 
Hugoniot  relations  to  appropriate  equation-of- 
state  data  for  sea  water  in  order  to  provide  (a) 
tables  of  the  desired  functions  up  to  very  high 
pressures  («*.  80  ldlobars)  for  use  in  the  theory 
of  propagation  of  underwater  explosion  waves,1 
and  (b)  tables  of  particle  and  propagation 
velocity  at  fairly  close  pressure  intervals  in  a 
lower  pressure  region  (up  to  ca.  14  ldlobars). 

n.  outuhs  ov  tbx  thiort  arc  com¬ 
putational  PROCXD0KXS 

tn  this  section  we  give  an  account  of  the 
hydrodynamicai  and  thermodynamical  relations, 
and  the  computational  procedures  leading  to  the 
numerical  remits  tabulated  in  Sections  III  and 
IV.  For  the  convenience  of  the  reader  a  glossary 
of  symbols  is  presented  in  Appendix  111. 

When  a  shock  wave  advances  with  velocity  U 
into  a  stationary  fluid  of  unperturbed  pressure 

*J.  C.  Kirkwood  and  H.  Bethe,  Th»  Prtsnrt  Win 
Produced  by  an  UndrmoUr  Expiation  ( Dept,  of  Commerce 
Bibliography  No.  PB  32182),  OSRD  Report  No.  388, 
Pert  1. 


S  .  IN?  H  -  .IsON 

Pi  and  specinc  voiume  vt.  the  pressure  p.  specific 
volume  ?,  and  particle  velocity  u  of  the  fluid 
behind  the  shock  front  are  determined  by  the 
Ran kine1- Hugoniot4  conditions,  which  express 
the  conservation  oi  mass,  momentum,  and 
energy  of  an  element  of  fluid  passing  through  the 
front.  For  the  purposes  of  this  paper,  these  con¬ 
ditions  may  conveniently  be  written 

u~Zip-pt)<vt-v)y,  .2.1) 

V  •t’sC(p-p»i/(®«-»)2*.  (2.2) 

.iff  -  (  1  /  2)  (p  -  pi)  (v-J-Vi) .  1 2.3) 

In  the  last  equation,  Iff  is  the  specific  enthalpy 
increment  of  an  element  of  fluid  when  it  passes 
through  the  front.  The  specific  enthalpy  is  de¬ 
fined  as  the  sum  of  the  internal  energy  per  gram 
and  the  pressure-volume  product,  pe. 

Given  equation  of  state  and  specific  heat  data 
for  the  fluid,  any  three  of  the  variables  p.  v,  C 
and  u  may  be  determined  as  functions  of  the 
fourth.  Here  we  shall  regard  p  as  the  independent 
variable.  For  certain  hydrodynamic  applications 
we  must  have,  in  addition  to  v,  U,  and  u  as 
functions  of  p.  the  sound  velocity 

c  — idp/dp)*j;  p— 1,®,  i2.4) 

the  Riemann  v-functiont 

f  IXp'.  Sj/cip',  S]]dp\  '2.5) 

•  »• 

and  the  undissipated  enthaipy 

—  f%  [p'.5]dp',  1.2.6) 

where  J  is  the  entropy. 

In  practice,  one  must  resort  to  successive 
approximations  to  effect  a  reduction  of  the 
Hugoniot  conditions,  combined  with  equation-of- 
state  and  specinc  heat  data,  to  a  set  of  relations 
expressing  u,  U,  and  r  as  functions  of  p.  To  this 

1 W.  J.  M.  Rankine.  Trans.  Rov.  Soc.  London.  AIM. 
277  f  1870). 

*  H.  Hugoniot,  J.  da  t'ecoie  poivt.  3t,  3  .2887);  S8.  I 
'1888). 

t  The  Riemann  »•(' unction  occurs  in  Rienunn’s  form  of 
the  hydrodynamicai  equations,  which,  for  the  case  ot 
jpherital  symmetry,  may  be  written  see  reference  U: 

(>  /SI)  +  it ■+• » )  H  lit  Hfr +•» » “  —  2 cu/r, 
fS/SI)  —  ■£  —  »)iS, Sri  j(»— a)— 0, 

where  l  is  the  time  and  -  is  che  radial  coordinate.  The 
other  quantities  have  already  been  defined. 


v. 
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end,  it  is  expedient  first  to  consider  certain 
quantities  as  functions  of  pressure  and  tem¬ 
perature.  ?  and  T.  or  pressure  and  entropy,  p  and 
5.  Before  proceeding  to  a  more  detailed  dis- 
cussion  of  the  calculations  it  may  perhaps  help 
to  orient  the  reader  if  we  consider,  qualitatively, 
contours  of  some  pertinent  quantities  in  the  p-T 
plane.  In  Fig.  1.  the  possible  states  of  a  given 
fluid  just  behind  the  shock  front  lie  along  a  single 
curve,  which  we  have  labeled  “Hugoniot."  An 
element  of  fluid  initially  in  the  state  [pi.  2  „) 
which  has  attained  a  state  {p,  T)  just  behind 
the  shock  front  finally  returns  to  a  state  {pi.  T {) 
along  the  adiabatic,  so  labeled  in  'h  »  figure.  Also 
included  are  the  designations  of  a  few  points  on 
a  p-S  basis  using  square  brackets  .v  'O'ding  to 
the  convention  introduced  a  Part  of  this 
section.  In  general.  T t  is  largar  than  T »  because 
of  the  dissipation  occuring  at  the  front.  The 
central  part  of  our  problem  is  the  determination 
of  the  Hugoniot  curve. 

We  shall  consider  in  Part  a  the  calculations 
due  to  Arons  and  Halverson4  which  are  intended 
to  be  accurate  in  the  range  of  relatively  low 
pressure  (ca.  0  to  20  kiiobars).  These  results  as 
stated  in  the  introduction  are  intended  for  the 
determination  of  the  peak  pressure  of  a  shock 
wave  from  measured  values  of  the  shock  front 
velocity  V  or  particle  velocity  u.  In  Part  b,  we 
shall  consider  the  calculations  of  Kirkwood  and 
Richardson,*  the  results  of  which  were  originally 
intended  for  the  applications  of  the  shock  wave 
propagation  theory  of  Kirkwood  and  Be  the4 
which  required  data  over  a  higher  range  of  pres¬ 
sures  {ca.  20  to  50  kiiobars). 

a*  Calculations  of  Arons  and  Halverson 

Here  we  outline  the  calculations'  suitable  for 
the  relatively  low  pressure  range  (ca.  0  to  20 
kiiobars)  based  upon  the  equation-of-state  and 
specific  heat  data  discuseed  in  detail  in  Appendix 
I.  For  the  range  0  to  1.5  kiiobars,  the  Ekman 
equation -of -state  was  used ;  in  the  range  0  to  25 

•A.  B.  Arons  sod  R.  R.  HsWsrson.  Huremet  Calaii s- 
tiami  /or  5m  Wa Ur  at  Urn  Shock  Prod,  OSRD  Report 
No.  6377,  NDRC  No.  A-469. 

*J.  G.  Kirkwood  and  J.  M.  Richardson,  Tho  Proaton 
Worn  Product*  by  an  Underwater  Exptenon.  Part  III.  OSRD 
Rsport  No.  *13  (Dept,  of  Commerce  Bibliography  No.  PB 
32114). 

’J.  G.  Kirkwood  and  E.  Mon  troll,  Prtmre  Warn 
Produced  by  an  Underwater  Expiotten,  II.  OSRD  Report 
No.  674  (Dept,  of  Commerce  Bibliography  PB-321S3). 


kiiobars.  the  Tait  equation-of-state. 

--•(O.  T)—  t’ip,  D  —  (l/n)  log£t-rp.'5f*)j, 

<-<r-273.16)8C.  (2.7) 

In  the  first  case,  the  initial  temperature  was 
.'a*  15°C;  in  the  second,  25°C.  In  both  cases, 
the  initial  pressure  po—0.  Neither  of  the  two 
equations-of-state  are  complete  in  the  sense  that 
T)  must  be  determined  by  auxiliary 
thermal  expansion  data  also  discussed  in  Ap¬ 
pendix  1). 

We  express  the  enthalpy  and  volume  incre¬ 
ments 

±H-H(p,  T)  ~H{ 0.  T,).  2  8) 

lvu>v(p,  T)—  pi.0.  To), 

in  terms  of  line  integrals,  first  along  an  isobar 
from  (0,  Ti)  to  ,0.  T)  and,  secondly,  along  an 
isotherm  from  (0,  7 *)  to  ip,  T)  (see  Fig.  1).  For 
the  enthalpy  increment  we  obtain 

Aff-A,tf-rArtf. 


c,(o.  r)dr-c^r. 


Artf- 


r  n- 


mp'.  n 
r . - 

ST 


t2.9) 


a  r-r-r,. 

where  c,( 0,  T)  is  the  specific  heat  extrapolated  to 
aero  pressure  and  tr  is  the  mean  of  c+  over  the 
temperature  range  A  T. 

For  the  volume  increment  we  obtain 

A**  A^fi-Arf, 

A^-e(0,  D -v(0,  r.)  -AsAT.  (2.10) 
Af»-etp,  T)— v(0,  T), 

where  $»  is  the  mean  thermal  expansion  at  zero 
pressure  over  the  temperature  range  A T. 

From  the  last  Hugoniot  condition.  Eq.  2.3), 
and  Eq.  (2.9)  we  obtain 

0(0,  r.)+(l/2)(A^)>-ArH 
Al*-— - - - ; - .  \.2.l l) 

cp-(i/2)(to) 

where  Ar£f  is  to  be  calculated  by  means  of  the 
third  of  Eq.  (2.9)  and  the  appropriate  equadon- 
of-state,  and  where  Ay*  is  to  be  obtained  from 
compressibility  data.  The  right-hand  side  of  Eq. 
(2.11)  depends,  of  course,  on  the  temperature  T. 
The  determination  of  A  T  is  accomplished  by  the 
method  of  successive  approximations  A  trial 
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Table  I.  N— ,  <  i.  md  iH  In  low  pressure  'ecion 
ba«fa  m  Efcman  -sjuacion-oi-itace  .  Sea  water :  fnitiai 
-.emperature  1 5 =  C  ;  <ai»nitv  52  parts  per  thousand  3.7d 
wt.  percent  NaCI);  :ii«4922.3  ft  -ecw  1500.3  m  iec. 


a  a 
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tctwr< 
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CO 

A 

m. 
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Ui 
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00 
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2 JO 
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12.7 

1.25 
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74.4 
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10.000 

5.48 
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15.0 

1  JO 

1 1.44 

57. 7 

34.1 

12.000 
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166.0 

10  0 

.  14.000 

744 

101.4 

22.1 

<  16.000 

544 

217.1 

23.2 

:  ts.ooo 

040 

242.0 

284 

,  20.000 

1047 

2664 

31.4 

{  22.000 

1143 

200.6 

34. J 

Table  11.  C—r.,  c«,  a.  and  A//  in  intermediate  pressure 
region  l.o  to  25  kilo  bars)  based  on  Tait  equation-ot. 
*tate:  **>7.900,  3~3.012).  Sea  water:  Initial  temperature 
25*C;  salinity  32  parts  per  thouand  1.3.79  wt.  percent 
VaClh  5014. 7  ft/ sec  «  1328.3  m/sec 
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1 
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J24 
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3340 
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UO.l 

50.000 

36.0 
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6.0 

3840 

277  3 

131.2 

! 

TO  000 

39.4 
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133.0 

5.0 

48.13 

346.1 

172.9 
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100.000 

424 

1014 

150.0 
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56.73 

406.0 

214.0 
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120.000 
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12.0 

64.70 

466.0 

234.8 

i 

140.000 

55.2 

1312 
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1 4.0 

72.25 
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203.2 
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•- 40.000 
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value  of  AT  i*  used  in  evaluating  the  right-hand 
side  giving  a  more  accurate  value  of  AT  on  the 
left-hand  side,  and  the  proceae  is  repeated  until 
the  results  of  two  successive  steps  differ  by  a 
sufficiently  small  amount.  One  or  two  steps 
generally  suffice. 

We  have  thus  obtained  T  as  a  function  of  p 
along  the  Hugoniot  curve  (see  Fig.  1).  It  is  now 
possible  to  calculate  immediately  the  particle 
velocity  u,  the  propagation  velocity  C,  and  die 
specific  volume  v  as  functions  of  p  behind  the 
shock  front.  The  results  using  the  Ekman 
equation-of-state  and  the  Tait  equation-of-state 
are  tabulated  in  Tablet  I  and  II,  respectively,  of 
Section  III. 

b.  The  Calculations  of  Richardson  and  Kirkwood 

Here  we  outline  the  calculation**  intended  for 
the  applications  of  the  shock  wave  propagation 


theory  oi  Kirkwood  and  Berne. ;  These  are  basee 
upon  the  equation-of-state  and  specific  heat 
data  discussed  in  detail  in  Appendix  II.  We 
use  a  modified  Tait  equauon-of-state  connecting 
vip.  7")  and  r<0.  T)  to  be  discussed  below.  In 
most  respects,  the  data  is  made  to  fit  the  proper¬ 
ties  of  an  aqueous  0.7  motai  NaCI  solution 
assumed  to  be  roughly  equivalent  to  sea  water 
of  salinity  r  -  32  parts  per  thousand  see  Section 
1  of  Appendix  I). 

In  these  calculations  the  initial  pressure  p»  is 
taken  to  be  zero,  and  several  different  initial 
temperatures  Ts>  are  used:  0°C,  20*C.  and  409C. 

Before  indicating  the  precise  nature  of  the 
modification  of  the  Tait  equation,  it  is  desirable 
to  mention  that  in  this  part  two  different  pairs 
of  independent  variables  will  be  used:  pressure 
and  temperature  f p.  T),  and  pressure  and 
entropy  ~p.  ST  Consequently,  in  order  to 
indicate  which  pair  are  used  in  a  function,  we  will 
use  parenthesis  to  indicate  the  first  pair  and 
square  brackets  to  indicate  the  second,  i.e. 
v(p.  T)  and  v\j >,  S]. 

The  modified  form  of  Tait  equation  introduced 
by  Kirkwood*  *  is 

!ogi»i/»)»(.l,  n)  log(,l-rP'-/l[S]).  -  2.12) 

where 

vmv£p,  5]—  rip,  7~J>,  5]),  t’i  S', 

(see  Fig.  1)  n  is  an  empirical  constant,  and  the 
function  .4CS]  is  related  to  the  function  3(f)  in 
the  original  isothermal  form  of  the  Tait  equa¬ 
tion.  Eq.  (2.7),  as  follows, 

i4CS]-3«C0,5]),  f-(r-273.l6)*C,  (2.13) 

The  reasons  for  introducing  this  modification  of 
the  Tait  equation  are  at  least  twotoid:  (1)  the 
anomaly  of  a  vanishing  specific  volume  vip.  T)  at 
a  finite  pressure  along  a  given  adiabatic  which 
does  not  differ  markedly  from  the  Hugoniot  curve 
in  the  case  of  water)  is  removed  to  a  higher  pres¬ 
sure  by  replacing  [c(0,  T)  —  vip,  T)2,  [tKO,  T)j 


Table  III.  Values  of  U —c./c.  tor  different  temperature* 
and  salinities  at  a  shock  wave  pesk  pressure  oi  1.00 
kilofasr. 


Seiiaity 

Ton  panuin 

U  — c*/<# 

(pans  per  (0001 

(•Cl 

i«V. 

32 

IS 

7.31 

32 

25 

7.31 

35 

15 

7.76 

17 
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by  loswCO.  51  vf>.  SJ.  and  2)  the  calculation 
of  quantities  denned  by  line  integrals  along 
adiabatics  is  greatly  simplified  by  talcing  S 
instead  of  T  as  one  of  the  independent  variables. 

The  (unction  .4C-S1  is  related  simply  to  ci,  the 
sound  velocity  at  zero  pressure  and  entropy  5 
according  to  Eq.  l2.4)  as  follows 

ri-c[0.5].  <2.U) 

On  the  basis  of  Bridgman’s  p-e-T  data  for  pure 
water,  an  average  value  ot  n  equal  to  7.  IS  has 
been  selected  for  the  present  calculations.  In 
Section  2  of  Appendix  II,  it  is  shown  that  n 
deviates  from  this  value  by  less  than  4  percent 
in  a  large  pressure-temperature  held  bounded  by 
adiabatics  starting  at  zero  pressure  and  tem¬ 
peratures  of  20*C  and  60*C,  respectively,  and 
extending  to  pressures  of  25,000  Icg/cra*.  We 
assume  that  n  has  the  same  value  for  an  aqueous 
0.7  motal  NaCl  solution  as  for  pure  water,  and 
we  obtain  by  interpolation  the  required  values 
of  B(t)  from  R.  E.  Gibson's  values  of  B(t)  for 
dilute  aqueous  NaCl  solutions  (see  Appendix  II, 
Section  1).  The  appropriate  heat  capacity  and 
thermal  expansion  data  are  discussed  in  Section 
1  of  Appendix  II. 

We  now  proceed  to  the  calculation  of  the 
quantities  u,  U,  c,  <r,  and  u.  We  first  express  these 
quantities  with  use  of  Eq.  (2.12)  in  terms  of  p. 

v~v(p,  T)  ~v£p,  5],  r«), 

m-e(0,  TO  -t0),  5],  and  ct-c(0.  T i)  «<0,  J] 
(see  Eqs.  (2.1H2.4),  also  Fig.  1)  as  follows: 


« 

■ 

1 

(2.15) 

U-pu/*, 

(2.16) 

c-c,(s,  /»)<-•»*. 

(2.17) 

*  «— — C(st/s)  (•-»<*  —  1 1 

IS  — 1 

Cl’ 

(2.18) 

/p)»-l  -1], 

n  —  1 

(2.19) 

Ones  the  temperature  T i,  to  which  an  element 
of  fluid  returns  along  the  adiabatic  intersecting 
the  Hugonioc  curve  at  [p,  T),  is  determined,  all 
of  the  above  quantities  may  be  determined  as 
functions  of  p.  To  accomplish  this,  the  enthalpy 
increment,  A/T,  occurring  in  the  third  Hugoniot 
condition,  Eq.  (2.3),  is  written  as  the  sum  of  two 
line  integrals,  the  first  along  an  isobar  from 


0.  Tu  “[0,  Stj  to  0.  rj-[0.  s:  and  the 
second  along  an  adiabatic  from  >0,  T,i  —  *0.  S] 
to  ( p,  T'i-Tp.S'  1  see  Fig.  1),  giving: 

Aff 


o,  Ddr, 


2.20) 


where  w  is  the  undisaipated  enthalpy  already 
defined  by  Eq.  (2.6)  with  0  and  given  ex¬ 
plicitly  in  terms  of  ci,  v\,  and  v  in  Eq.  (2.19). 
The  dissipated  enthalpy  k  can  be  determined  as 
an  explicit  function  of  T»  and  T i  from  specific 
heat  data  (Appendix  II,  Section  1).  Combining  the 
third  Hugoniot  condition.  Eq.  (2.3),  with  Eqs. 
(2.19)  and  (2.20)  we  obtain  the  relation 


- O-l).  (2.2) 

2nt>i 


Tasls  IV.  Properties  of  tea  water  at  a  tirade  front. 
(Initial  temperature  0*C;  salinity  0.7  m  NaCl;  C»— 1+43 
m/tec.) 


P 

» 

u 

e 

# 

*X!0-« 

h 

f 

UJo» 

(»/ 

(m/ 

m/ 

(m/ 

(m  / 

(joule/ 

(cm*/ 

twr> 

ne) 

MC) 

MC) 

MC) 

HCf* 

cm) 

091 

0 

0 

__ 

0 

0 

0 

0.9915 

5 

257 A 

19J0 

2140 

255.5 

0.4560 

6.740 

*395 

10 

455A 

2399 

2720 

4204 

OA720 

2540 

.0040 

IS 

S7SJ9 

2505 

3145 

352A 

1.270 

34.40 

.7710 

20 

ffii 

2045 

3510 

664A 

1-655 

5645 

.740J 

25 

5075 

5055 

7694 

2-050 

1224 

.7519 

JO 

90S  A 

5205 

4125 

5554 

2.405 

1604 

.7106 

is 

99TA 

3400 

4J4S 

9404 

2.770 

2014 

.791$ 

40 

1000 

36M 

4640 

1020 

3.140 

244A 

.6909 

SO 

1240 

toes 

3095 

1175 

3440 

351.0 

.6042 

60 

1511 

4500 

5465 

1515 

♦475 

419A 

.6720 

70 

1515 

4505 

SOTO 

1455 

5-2*5 

5090 

.6641 

O 

tus 

455 

6225 

1505 

4.000 

3954 

.6579 

90 

1740 

5120 

6570 

1705 

6.750 

6764 

6542 

Tasls  V. 

Proper  ties  of 

6061  water  ot 

a  shock  front. 

(Initial  temperature  20  *C ; 

salinity  0.7  m  NaCl;  C# 

-1317 

m/see.) 

(kSw 

• 

V 

e 

9 

*XlO-» 

k 

9 

< ■ < 

(t*/ 

(*/ 

(*/ 

[mt 

(jo UM/ 

(cm*/ 

b m) 

Me) 

Mel 

MC) 

MC)* 

09) 

f«*) 

0 

0 

0 

0 

0 

0.9929 

.  s 

23IA 

197 1 

22JO 

2404 

0.4595 

5470 

4660 

to 

4254 

nu 

2755 

4154 

04790 

25.45 

4120 

15 

567  A 

24JO 

3175 

S44.0 

1400 

♦943 

.7707 

20 

099  JO 

2000 

3555 

665.0 

1.670 

SO  AS 

.7555 

IS 

740.0 

3110 

3055 

76SA 

2.050 

USA 

.7501 

JO 

0974 

3520 

4140 

554.0 

2.425 

2.7VS 

1524 

.7245 

35 

940.0 

3510 

4405 

9464 

1V2A 

.7150 

40 

t07S 

3000 

ms® 

1050 

3.160 

255A 

7054 

30 

1255 

4030 

5100 

1105 

1405 

3174 

.6600 

60 

1500 

4325 

ISOS 

1550 

4.605 

♦044 

.6765 

70 

1510 

4610 

3000 

1465 

SJtO 

♦094 

6679 

10 

1625 

tats 

6260 

1600 

6A50 

3734 

.6626 

18 
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Tails  VI.  Properties  ol  sea  water  at  a  shock  front. 
Initial  temperature  +0°C,  salinity  0.7  m  NaC!.. 


9  a  li  i  •  *X  10**  k  f 


lulo- 

m/ 

m/ 

•  m/ 

m- 

cm*' 

tMuri 

me) 

wc; 

me) 

mc) 

me  i» 

(mi 

zm> 

0 

9 

_ 

J 

0 

0 

0.9993 

3 

249.5 

2005 

2253 

247.0 

0.4410 

549 

*749 

10 

421.5 

2140 

277S 

415.0 

04*70 

22.75 

319* 

SI 

5664 

IMS 

1193 

550.0 

1.290 

4*40 

7*59 

20 

6*9.0 

2900 

3550 

ViAi 

1  o*5 

’*.70 

.762! 

25 

7964 

3110 

JM5 

7704 

2.065 

!  134 

7441 

.10 

599.0 

1135 

41JO 

346.0 

2.445 

15  to 

.729* 

15 

992.0 

1525 

3705 

4415 

9534 

1040 

2415 

1394 

7179 

40 

tooo 

4440 

3.1*5 

2104 

70*0 

SO 

1240 

4035 

5110 

1200 

3.915 

J15.0 

4926 

00 

1140 

4140 

5515 

IMS 

4.440 

400.5 

4*11 

TO 

tsio 

4415 

1900 

1455 

5470 

4*1.0 

.6717 

where  y  » (’»»,’»)  *.  With  the  aid  of  tables  of  a/ci5 
and  oi  as  functions  of  7\  and  Tt,  Eq.  '.2.21)  may 
be  solved  by  successive  approximations  giving  T i 
as  a  function  of  the  parameter  y.  Since  the  equa¬ 
tion  of  state,  Eq.  (2.12),  may  be  expressed 
simply  as  pmB(T:  —  273.16)(jr—  lj,  the  tem¬ 
perature  T i  may  be  determined  as  a  function  of 
the  pressure  p  by  a  tabular  elimination  of  y.  By 
graphical  interpolation,  T i  is  finally  determined 
for  the  desired  integral  values  of  p  (in  kiloban), 
and  the  functions  a,  U,  c,  <r,  and  u  are  then  com¬ 
puted  as  functions  of  p  by  means  of  Eqs.  (2.15)- 
(2.19). 

OL  NUMERICAL  RESULTS  OR  ARONS  AND 
HALVERSON 

in  fundamental  shock  wave  studies,  it  is  fre¬ 
quently  necessary  to  know  values  of  U —c%/c% 
and  u  at  given  pressure  levels  to  the  highest 
possible  degree  of  accuracy.  With  this  object  in 
view,  the  calculation  methods  described  in 
Section  2A  were  applied  to  the  best  available 
equadon-to-state  data.  The  numerical  results  are 
given  in  Tables  1  and  II.  A  critical  discussion  of 
the  equation-of-state  data  will  be  found  in  Ap¬ 
pendix  I  together  with  references  to  the  sources 
from  which  they  were  obtained- 

Table  I  gives  results  for  the  “low  pressure" 
region,  covering  shock  wave  peak  pressures  of 
from  0  to  1.50  kiloban  (ca.  22,000  p.s.i.).  The 
calculations  in  this  table  were  baaed  upon  the 
Ekman  equation -of -state  for  sea  water  (see  Ap¬ 
pendix  I)  which  is  used  in  the  calculation  of 
sound  velocity  for  echo-ranging  tables. 

Since  the  Ekman  equation  deviates  appreci¬ 
ably  from  experimental  compressibility  data  at 
pressures  exceeding  2  kilobars.  this  equation  was 
abandoned  in  the  “intermediate  pressure"  region. 


The  results  in  Table  1 1  are  applicable  principally 
to  the  region  between  1.5  and  14  kilobars  ca. 
200,000  p.s.i.)  and  are  based  on  a  careful  fit  of 
the  Tait  equation  to  Adams's  experimental  com¬ 
pressibility  data  'see  Appendix  1). 

Tables  1  and  II  were  computed  for  certain 
specific  values  of  temperature  and  sea  water 
salinity  (equivalent  to  0.675  moial  N’aCl),  and 
it  is  shown  In  Table  111  that  the  value  of 
U —Ci/  Ct  is  not  very  sensitive  to  changes  in  these 
variables. 

IV.  NUMERICAL  RESULTS  OR  K2JUCWOOD 
AND  RICHARDSON 

In  Tables  IV'  to  VI,  the  particle  velocity  u,  the 
shock  front  velocity  U.  the  sound  velocity  c,  the 
Riemann  »-f unction,  the  undissipated  enthalpy 
w.  the  dissipated  enthalpy  k ,  and  the  specific 
volume  e  of  sea  water  (0.7  molal  NaCl  solution) 
are  presented  as  functions  of  pressure  p  along 
three  Hugoniot  curves,  starting  at  zero  pressure 
and  the  temperatures  0°C»  20°C  and  40“C. 
respectively.  These  results  have  been  calculated 
by  the  procedures  of  Part  b  of  Section  2  and 
the  data  of  Appendix  II.  The  results  above  30 
kilobars  represent  extrapolations  beyond  the 
range  of  experimental  data:  consequently  the 
validity  of  the  results  above,  say,  50  kilobars,  is 
questionable. 

In  closing  this  discussion  of  the  calculations, 
the  authors  wish  to  acknowledge  their  gratitude 
and  appreciation  to  Professor  J.  G.  Kirkwood  of 
Cornell  University  for  his  contributions  in 
initiating  the  work  and  in  supplying  valuable 
guidance  and  advice. 

APPENDS  Iff 

1.  Salinity  and  Temperature  Conditions 

All  calculations  were  made  for  sea  water 
having  a  salinity  of  32  parts  per  thousand  vthe 
average  salinity  it  sea  water  at  Woods  Hole, 
Massachusetts).  Salinity  is  defined  in  terms  of 
directly  measured  chiorinity  as : 

s -0.050+ 1.8050  Cl 

where  s  and  Cl  are  expressed  in  parts  per 
thousand. 

It  was  calculated  from  the  average  composition 

tt  Equation -of -tut*  data  uaed  in  computation  of  • 
Tablet I  and  II. 
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of  sea  water  that  .on  the  basis  of  ionic  strength)  a 
salinity  of  32  parts  per  thousand  is  equivalent 
to  an  NaCl  solution  having  a  molality  of  0.675 
or  a  weight  percentage  of  3.79  percent  NaCl. 

Table  1  was  computed  for  an  initial  tempera* 
cure  of  15*C  because  this  temperature  is  a  rough 
average  of  conditions  normally  encountered  in 
experimental  work.  Table  II  was  computed  for 
an  initial  temperature  of  25*C  because  this  was 
the  temperature  quoted  for  the  available  com- 
presaibility  data.1  Table  III  shows  that  the 
results  are  not  sensitive  to  small  variations  in 
temperature  and  salinity. 

2.  Specific  Volume  and  Coefficient  and 
Thermal  Expansion 

The  best  sources  of  data  seem  to  be  the 
oceanographical  rabies  of  Knudsen.*  Second 
power  equations  in  t(*C)  were  fitted  to  the  data 
tabulated  for  s  -  32 : 

For  Table  I : 

s(/)  —  0.97709 + 2.05  X 10 — *(t  — 15) 

+4XlQ“4(/— 15)*. 

For  Table  II: 

v(t) -0.97956+2.85  X  10“*(<-  25) 

+4Xl0-»(t-25)». 

3.  Heat  Capacity 

The  beat  capacity  data  used  in  computing 
Tables  I  and  II  are  those  quoted  by  S.  Kuwa- 
hara:1* 

C9-C„ -0.0004226/+0.000006321/*  caL/gm'C. 

Tasls  VII.  Comparison  of  eipsrimsn  tally  measured 
sound  velocity  with  calculations  baasd  on  tho  Ekraan 
Com  pro— I  hili  ty  Equation.  (Salinity— 31.7  parts  per 
thousand  ) 


T'T-ST1* 

Velocity  of  «*ad  (ft/*ec> 
Mworori  Cokiiiwod 

Drristioa 

(T»l 

10.9 

4887.7 

4673.6 

0.26 

11.6 

1UU 

4883.8 

0.04 

U.6 

4893, 1 

4883.8 

0.19 

11.5 

4902.4 

48813 

0.41 

11.1 

4MU 

4878.5 

OJO 

•  Adams.  J.  Am.  Chen.  Soc  SI,  3769  (1931). 

'  Qcsassrephifsi  Tobin,  Comimariat  o f  Agriculture, 
USSR,  Moscow,  1931.  (A  general  compilation  of  ocsnno* 
graphic  data  by  N.  N.  Zubov.) 

“  S.  Kuwahara,  Vdomtj  of  Sound  in  Sou  Wator  and 
Caleuiatiom  of  tho  Vtionty  for  Us*  in  Some  Saundtn* 
(Hydrographic  Dept.  I.J.N.  Tokyo,  193». 


Table  VUI.  Comparison  oi  \damst  experimental 
compressibilities  and  the  empirical  he  given  by  the  Elcman 
and  Tait  equations. 


P 

■  kber  > 

Adorn  » 
experimental) 

Pure  3.79% 

HO  VeC T 

Ekata 

equation 
S-32 
(T able  n 

Ta*t«  equation 

8  — J.Q12  h~7  44$ 
Table  ID  S-J.Ud 

0.0 

0.0000 

0.0000 

0.0000 

o.oooo 

0.0000 

0.3 

.0212 

.0196 

.0198 

■  .0197 

0198 

!  0 

.0393 

.0368 

.0370 

0368 

.0370 

1.3 

.0535 

.0322 

.0522 

0518 

.0522 

2.0 

i  0699 

.0658 

.0653 

.0633 

0659 

3.0 

0945 

.0894 

.0871 

.0887 

.0897 

4.0 

■  .1152 

.1091 

— 

.1083 

.1095 

5.0 

.1330 

.1263 

— 

.1234 

.1275 

6.0 

1  .1485 

.1417 

— 

.1405 

.1431 

7.0 

'  .1622 

.1532 

— 

.1540 

.1569 

3.0 

'•  .1746 

.1670 

— 

.1662 

.1695 

9.0 

!  .1858 

.1781 

— 

.1775 

.1812 

10.0 

1  .1964 

.1886 

— 

.1876 

.1917 

11.0 

j  .2059 

.1980 

“ “ 

.1972 

.2017 

*  l,«!  lo«<l  +S/S1. 


where 

Cw  - 1.005  — 0.004136i+0.0001098i* 

—0.0000013245* 

In  the  above  equations,  t  is  temperature  in 
and  s  is  salinity  in  para  per  1000-  These  data 
are  in  good  agreement  with  those  used  by  Kirk¬ 
wood  and  Richardson,  quote  in  Appendix  II. 

4.  Compressibility  Date  for  Low  Pressure  Region 
(Table  I) 

The  following  equation  was  used  in  computing 
Table  I : 

4886 

10»tt  - £227 + 28.33/  -0.551/* 

l+0.183p 

+0.004/']+p[10S.5+9.50/-0.158/*] 

(t  —  28\ 

-  1.5p*/- [147.3 -2.72/ 

+0.04/*  -  j>(32.4  -  0.87/ +-0.02/*)  ] 
/y-28\  5 

+  ^-^-j  [4.5— 0.1/— p(  1.8— 0.06/)], 

where  p  is  pressure  in  Idiobara,  t  is  temperature 
in  degrees  centigrade,  and  a  is  defined  by : 

r-pt(l -up), 

y  is  defined  by : 

7- -0.069+1.4708  Cl— 0.001570  Cl* 

+0.0000398  Cl*. 
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The  above  empirical  equation  for  sea  water 
compressibility  is  due  to  Ekman-1  and  has  been 
widely  used  for  computation  of  sound  velocity 
in  sea  water. 14  u  The  validity  of  the  Ekman 
equation  for  sound  velocity  calculations  was 
verified  experimentally  as  indicated  in  Table  VII. 

Experimental  sound  velocity  measurements 
were  made  by  recording  with  a  rotating  drum 
camera  the  signals  applied  to  a  cathode  ray 
oscilloscope  by  two  very  small  piezoelectric 
gauges  placed  a  known  distance  apart.  The  sound 
source  was  -  No.  8  detonator  cap  placed  far 
enough  away  from  the  gauges  so  that  the  effect 
of  finite  pressure  amplitude  was  less  than  0.03 
percent.  An  error  of  about  0.2  percent  was  in¬ 
herent  in  the  experimental  work  owing  to  slight 
errors  in  the  alignment  of  the  two  recording 
gauges  with  the  sound  source.  This  accounts  for 
the  magnitude  and  systematic  nature  of  the  dis¬ 
crepancy  apparent  in  Table  VII. 

Further  verification  of  the  applicability  of 
Ekman’s  equation  in  the  region  up  to  1.30 
kilobars  is  given  in  Table  VIII  where  values  ob¬ 
tained  from  the  equation  are  compared  with  the 
experimental  values  of  Adams  for  NaCl  solutions. 

5.  Compressibility  Data  for  Intermediate  Pres¬ 
sure  Region  Table  II) 

As  indicated  in  Part  1  of  this  appendix,  a  sea 
water  salinity  of  32  parts  per  thousand  cor¬ 
responds  to  a  3.79  weight  percent  solution  of 
N'aCl.  The  compressibility  of  NaCl  solution  of 
this  concentration  was  obtained  by  graphical 
interpolation  of  Adams's  data.* 

The  Tait  equation  in  the  form  ■ 

®(0,  T)-v(p,  DM 0.  r)-(l/n)  logCH-p/B(<)3, 
<-(r-. 273.16)#C 


was  then  fitted  to  Adams's  data.  In  an  effort  to 

Table  IX.  Value*  of  *  computed  from  p—v—T  data 
.ueing  n«  7.15  in  computation  of  XT'.. 


*• 

l«C) 

soot 

t  (*/«■*> 
15.000 

23.000 

20 

7.211 

7.13} 

7.150 

40 

7.5«0 

7.126 

6.969 

60 

7.411 

7.084 

6.86 8 

^£<MftMiam4e  ueovuut  ibraifk  *McH  tbe  a*i**atic  for  j 


UV.  W.  Eicman,  Pubiications  is  Circons  tones  So.  4 3 
(Conceil  Peimanent  Internationale  Pour  L'Exploratioo  de 
la  Mer.  November  1908). 

“  Matthew*.  Tabiss  of  tks  VtiocUj  of  Sound  in  Pun 
Worn  and  Son  Watsr  for  Us*  tn  Echo  Soundtnt  and  Sound 
Panftnf  (Hydrographic  Dept.,  Admiralty,  H.D.  No.  282). 


check  the  T.iit  equation  against  the  Ekman 
equation  used  for  computation  of  Table  I.  a  fit 
was  first  made  to  the  lower  pressure  region. 
Values  of  i*  and  Bi25°C)  were  so  selected  that  the 
equation  not  only  fitted  the  data  of  Adams  with 
adequate  precision  but  also  yielded  the  correct 
velocity’  of  sound  in  the  limit  of  zero  pressure. 
This  additional  restriction  that  the  equation 
give  c»  — 1528  m  sec.  at  2S°C  and  a  — 32)  required 
that  i»fJ125°C)  -  23.497,  the  latter  relation  being 
obtained  from  the  thermodynamic  equations: 


(dv  \  tv  T  /  d tr 
dp/  r  Co1  c?\  dT 

(do  \  V» 

—  }  - - at  p-0. 

dpJ  r  nB{t) 

In  this  case  *  was  taken  as  7.445  and  B(25’C) 
as  3.156  kilobars,  and  the  resulting  equation  fits 
the  data  of  Adams  quite  closely  up  to  pressures 
of  about  4  kilobars  as  shown  in  Table  VII  l.  For 
purposes  of  further  calculation,  the  temperature 
variation  of  B  was  assumed  to  be  the  same  as  that 
used  by  Kirkwood  and  Richardson  on  the  basis 
of  a  private  communication  from  Gibson  see 
Appendix  II).  Calculation  of  U—c<>,cq  at  1.00 
kilobar  yielded  a  value  of  7.85  percent,  in  good 
agreement  with  the  value  of  7.81  percent  ob¬ 
tained  from  the  Ekman  equation. 

Having  verified  the  accuracy  of  results  ob¬ 
tained  from  the  Tait  equation  when  fitted  as 
described  above,  the  same  technique  was  used  to 
fit  the  equation  to  the  intermediate  pressure 
range  (up  to -values  for  il  kilobars  quoted  by 
Adams;  it  was  assumed  safe  to  extrapolate  the 
resulting  equation  to  pressures  of  14  or  15  kilo¬ 
bars).  it  was  found  that  the  best  fit  of  the  data 
as  well  as  a  correct  value  for  the  velocity  of 
sound  were  obtained  by  taking  it  —  7.300  and 
B(25°C) -3.012,  the  temperature  variatior  of  3 
again  being  assumed  to  be  that  mentioned  above. 
The  Tait  equation  containing  these  parameters 
was  then  used  for  the  computation  of  Table  II. 
The  fit  of  the  equation  to  Adams's  data  is  shown 
in  Table  VIII. 
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APPENDIX  a 

1.  Data  Employed  ia  the  Computations  of  Part 
B  of  Section  H 

In  the  modified  Tait  equation,  Eq.  v2.12),  the 
function  -4[S]-B(<),  where  I-7T0,  53  —  273.16, 
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is  determined  from  the  empirical  values  of  S'O, 
titting  the  original  isothermal  Tair  equation,  Eq. 
>2.7),  to  experimental  data.  R.  E.  Gibson11  jives 
third-degree  .'-expansions  of  Bit)  for  various 
molalities  of  NaCl.  By  interpolating  the  coef¬ 
ficients  '.the  constant  term  numerically  and  the 
other  graphically)  for  a  molality  of  0.7.  one 
obtains  B<0  »3.134— 1-65  X  10~*(<  —  55)  — 1.181 
X  10“*(l  —  55),-ta5.32X  I0~r,l  —  55)*  Icilobars. 

The  specific  heat  crK 0,  D  for  a  0.7  molal  NaCl 
solution  was  obtained  by  Interpolation  from  the 
values  quoted  in  the  International  Critical  Tables 
and  Physikalischchemische  Tabeilen.  The  resulting 
set  of  values  is  fitted  adequately  by  the  ex¬ 
pression 

c,l0.  t-i- 273.16) -3.9644-p-6.24 

X  10~4t  joule/  gm.  deg. 

From  Gibson  and  Loetfler14  a  set  of  values  of 
v(0.  T)  covering  the  range  from  2 S*C  to  95°C 
inclusive  was  obtained  for  a  0.7  molal  NaCl 
solution  by  means  of  empirical  equations  giving 
riO,  D  as  a  function  of  concentration  for  each 
temperature.  In  extrapolating  to  higher  tem¬ 
peratures.  the  relation, 

s(0.  1+273.16)  -0.994150-4- 2.929  X  lO”4^  —  25) 

+ 3. 241 X  —  25) *  cm* / gm . 

was  used  ;  for  lower  temperatures  <  10°C), 

v(0,  1+273.16)  -0.991442+6.025 

X  10"*(t  -  3.8)1  em’/gm. 

Z.  Test  of  the  Modified  Tait  Equation  with 
Bridgman’s  Data  for  Pure  Water.  Deter¬ 
mination  of  the  Characteristic 
Constant  n. 

The  modified  Tait  equation-oi-state,  Eq. 
(2.12),  may  for  our  present  purposes  be  written 
in  the  form 

log(s0>,  5]/tO-  ■SD 

-(l/»)  Iog(l+p//4CJ7)  (II-l) 

where  is  related  to  the  B(t)  in  the  original 
isothermal  equation  of  state  as  follows, 

ACS]-3(lt), 

1*—  7"»  — 273.16,  (11*2) 

_  r,-ito.  5]. 

*  Priests  communication. 

u  Gibson  and  Lowlier,  J.  Am.  Cbwti.  Soc.  51.  443  >,l9*t). 
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According  to  the  convention  introduced  in  Part 
B  of  Section  II.  parentheses  '  )  after  a  function 
denotes  that  the  independent  variables  are  p  and 
r,  whereas  square  brackets  7  ’  denote  that  they 
are  p  and  S. 

Now  we  wish  to  test  Eq.  II-l)  with  Bridg¬ 
mans11  p-v-T  data  for  pure  water  with  the 
ultimate  object  of  linding  the  best  value  for  «. 
We  assume  implicitly  that  n  does  not  vary 
rapidly  with  NaCl  concentration.  To  make  the 
comparison,  we  first  must  know  the  values  of  the 
temperature  T  corresponding  to  the  various 
points  Z?~  'he  calculation  of  which  we  con¬ 
sider  below. 

Letting  F'0,  53—7%,  T*j>,  SZ  —  T,  and  T—T, 
mAT,  we  have 


ST- 


■I. 


■»  dTZp.  5] 

dp 


'■dp- 


I, 


»  dv[p,  5] 

IF 


dp.  ,11 


Using  Eq.  1 1 1- 11.  a  simple  calculation  yields 
G 

ST - - [(l-D)(l+#>  .4) 

a-b'.4)';’4 

II— 4) 

where 


A  —  <4(LSj  —  BUtu, 

.4'Cs>ro.5]  r»3'u,)-n o. 

Qm - — - , 

p»— i  r,) 


psa[5} 


Dm- 


dv[0.  Sj 
3S 


.4T5X0.  5] 


*5(1»)0» 

•B’Us/tpfO,  T i) 


r-ro 


To  calculate  T,  given  a  specified  p  and 
r»-7t0,  SI.  a  tentative  value  of  *—7.15  was 
chosen  for  use  in  Eq.  (11-4).  The  corresponding 
value  of  r(j».  5]-r(p.  21  was  obtained  by  inter¬ 
polation  from  Bridgman’s1*  p-v-T  data.  In¬ 
serting  these  values  of  vfp,  5]  in  Eq.  (II-l),  and 
knowing  the  values  of  r£0.  5}— tn,0,  To)  and 
3(U)  for  pure  water,  a  set  of  values  of  *  was 
calculated  for  p  —  5,000,  15.000.  25,000  kg,  cm* 

“  Bridgman.  J.  Chem.  Phv*.  J,  397  ;  1933)  md  private 
communication. 
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ind  :.)»  Tj  —  273.16*  20®,  40®.  oO®C.  me  resuits 

summarized  in  Table  IX. 

Additional  data  for  pure  water)  used  in  Eqs. 
11-1)  and  11-4)  were 

Bit )  -  2.996—7.285  X  -  25)  -  1.790 

X  If*-*'!  -  2S)*-t-6. 13  X  10~T'l  -  25)*  '<i!obars.“ 

and 

1  d  logri.O,  273.16— tl 
2.303 


2U-3.98) 

"  244.360+ 15, 040tt-3.98)«« 

l0.62)(15,040)(t-3.98)l“ 

"244,860 + 15 ,040<1  -  3 .98)*  **]‘ 

obtained  from  Ipatov's1*  empirical  equation  for 
v  by  differentiation. 

The  average  value  of  n  is  7.146.  In  the  present 
calculations  this  value  has  been  rounded  off  to 
7.15. 

The  entries  in  Tables  IV,  V,  and  VI  therefore 
contain  more  significant  figures  than  the  test 
justifies..  On  the  basis  of  the  test,  the  errors 
associated  with  the  use  of  the  modified  Tait 
equation  are  of  the  order  of  several  percent.  In 

«  I.  V.  Ipatov.  J.  Phyj.  Chetn.  U.S.S.ft.)  S.  1230  (1934). 


particular,  the  resuits  obtained  for  ow  pressures 
will  disagree  with  .mown  data  bv  several  percent. 

APPIKDIX  m 

Symbols 

i  [5’ —  parameter  in  modified  '"adiabatic  ' '  Taitequation- 
ii- » tat*. 

3U)  —  parameter  in  isothermal  Tait  equation-oi-scate. 

-  —  local  velocity  of  wound. 

veioary  of  wound  at  zero  pressure. 

—  s  peace  beat  at  constant  pressure. 

fc -dissipated  enthalpy  increment:  J ^  r,v0.  T'dV 

•enthalpy  increment:  .!#•<»— «. 

*  •characteristic  constant  in  Tait  equation-of -state, 
^•initial  pressure  ahead  of  shock  front.  p»— 0,  in 
these  calculations, 
p-  pressure  behind  shock  front. 

3-entropy, 
i  -sea  water  salinity, 
i  •  temperature  in  ’C. 

T -absolute  temperature, 
a  — particle  velocity  behind  shock  front, 
f— shock  front  propagation  velociry. 
ra— specific  volume  of  medium  ahead  of  shock  front, 
o  •  specific  volume  of  medium  behind  shock  front. 

3, •mean  compressibility  at  zero  pressure  over  temper- 
ature  range 
»•  density. 

9 -Riemann  function:  f  —r-—rzdp' 

•f"  c[fi'.  5] 

w  "L p-  djfp’. 
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appendix  a 

COMPUTER  LISTING 


01  L3LTHH 

26  STO  05 

51  X 

76  - 

02  TLOOP? 

27  ^PRESSURE? 

52  RCL  25 

77  STO  07 

03  PROMPT 

28  PROMPT 

53  - 

78  THH  • 

04  STO  34 

29  STO  12 

54  CHS 

79  ARCL  X 

05  TN? 

30  RCL  05 

55  RCL  25 

30  AVIEW 

06  PROMPT 

31  XEQ^VW 

56  RCL  23 

31  STOP 

07  STO  28 

32  STO  27 

57  Y*X 

32  RCL  06 

08  1 

33  USL  01 

58  - 

33  1000 

09  + 

34  XEQ^T 

59  2 

34  X 

10  RCL  28 

35  STO  25 

60  / 

35  THT  - 

11  1 

36  RCL  01 

61  RCL  20 

36  ARCL  X 

12  - 

37  XEQ^VW 

62  X 

87  AVIEW 

13  / 

38  STO  26 

63  RCL  26 

88  STOP 

14  STO  21 

39  RCL  01 

64  RCL  27 

39  RCL  07 

15  1 

40  XEQTSPEED 

65  - 

90  RCL  06 

16  RCL  28 

41  STO  20 

66  2 

91  1000 

17  1/X 

42  RCL  01 

67  / 

92  X 

18  — 

43  XEQTHT 

68  RCL  20 

93  - 

19  STO  22 

44  STO  06 

69  X 

94  THH-HT 

20  RCL  28 

45  RCL  25 

70  RCL  26 

95  ARCL  X 

21  1/X 

46  RCL  22 

71  / 

96  AVIEW 

22  CBS 

47  Y*X 

72  RCL  25 

97  STOP 

23  STO  23 

48  1 

73  1 

98  ABS 

T 

74  - 

99  RCL  34 

24  TEMP  07 

49  - 

25  PROMPT 

50  RCL  21 

75  X 

100  X<>Y 

25 


101  X<Y? 

102  GTO1 SHOCK 

103  GTO  01 

104  LBLTVW 

105  25 

106  - 

107  STO  03 

108  ENTER* 

109  X 

110  3.241E-6 

111  X 

112  RCL  03 

113  2.929E-4 

114  X 

115  + 

116  0.99415 

117  + 

118  RTN 

119  L3LTT 

120  TTEMP  1? 

121  PROMPT 

122  STO  01 

123  XEQTBT 

124  RCL  08 

125  1/X 


126  RCL  12 

151  ENTER* 

127  X 

152  X 

128  1 

153  RCL  05 

129  + 

154  ENTER+ 

130  STO  25 

155  X 

131  RTN 

156  - 

132  LBLTSPEED 

157  0.000312 

133  STO  01 

158  X 

134  XEQTBT 

159  + 

135  STO  11 

160  RTN 

136  RCL  01 

161  lbltbtprime 

137  XEQTVW 

162  ENTER* 

138  0.001 

163  55 

139  X 

164  - 

140  RCL  U 

165  STO  10 

141  X 

166  -2.362E-4 

142  STO  20 

167  X 

143  RTN 

168  -0.00165 

144  LBLTHT 

169  + 

145  STO  01 

170  RCL  10 

146  RCL  05 

171  ENTER* 

147  - 

172  X 

148  3.9644 

173  1.596E-6 

149  X 

174  X 

150  RCL  01 

175  + 

176  1.013E8 

201  X 

177  X 

202  RCL  02 

173  STO  24 

203  -1.13 IE -4 

179  RTN 

204  X 

180  LBLT8ETA 

205  + 

181  25 

206  RCL  04 

182  - 

207  -0.00165 

183  6.482E-6 

208  X 

184  X 

209  + 

185  2.929E-4 

210  3.134 

186  + 

211  + 

187  0.001 

212  STO  08 

188  X 

213  1.013E8 

189  STO  19 

214  X 

190  RTN 

215  STO  09 

191  LBLTBT 

216  RTN 

192  55 

217  LBLTSHOCK 

193  - 

218  RCL  11 

194  ENTER* 

219  1.013E8 

195  ENTER* 

220  / 

196  STO  04 

221  1/X 

197  X 

222  RCL  12 

198  STO  02 

223  X 

199  X 

224  1 

200  5.32E-7 

225  + 

26 


226  RCL  28 

251  TWATER  U  1 

227  1/X 

252  ARCL  X 

228  YtX 

253  AVIEW 

229  1/X 

254  STOP 

230  RCL  26 

255  RCL  15 

231  X 

256  RCL  27 

232  STO  13 

257  0.001 

233  FIX  6 

258  X 

234  h  - 

259  X 

233  ARCL  X 

260  RCL  14 

236  AVIEW 

261  / 

237  STOP 

262  STO  16 

238  FIX  0 

263  TSHQCX  0 

239  RCL  12 

264  ARCL  X 

240  1.013E8 

265  AVIEW 

241  X 

266  STOP 

242  STO  13 

267  RCL  20 

243  ACL  27 

268  RCL  28 

244  RCL  13 

269  X 

245  - 

270  SQRT 

246  0.001 

271  STO  17 

247  X 

272  RCL  26 

248  X 

273  RCL  13 

249  SQRT 

274  / 

250  STO  14 

275  RCL  28 

276  1 

301  X 

277  - 

302  RCL  27 

273  2 

303  0.001 

279  / 

304  X 

280  YtX 

305  X 

281  X 

306  RCL  28 

282  STO  18 

307  1 

283  ^WAVE  C  • 

308  - 

284  ARCL  X 

309  / 

285  AVIEW 

310  RCL  29 

286  STOP 

311  / 

287  GTOTOEL  T 

312  STO  30 

288  RTN 

313  RTN 

289  LBLTC? 

314  L3LTI)2 

290  0.624 

315  RCL  28 

291  X 

316  RCL  08 

292  3964.4 

317  X 

293  + 

318  RCL  19 

294  STO  29 

319  X 

295  RTN 

320  RCL  24 

296  LBLTG1 

321  1.013E8 

297  RCL  05 

322  / 

298  273.16 

323  / 

299  + 

324  RCL  27 

300  RCL  24 

325  0.001 

27 


326  X 


351  RCL  31 

352  + 

353  X 


327  / 


328  STO  31 

333  & 

329  RTS 

354  CHS 

330  LBL^TJEL  T 

355  RCL  28 

331  RCL  12 

356  + 

332  RCL  08 

357  1 

333  / 

358  - 

334  1 

359  RCL  32 

335  + 

360  RCL  31 

336  STO  32 

361  1 

337  RCL  05 

362  + 

338  XEQTCP 

363  X 

339  RCL  05 

364  + 

340  XEQTBTPRME 

365  RCL  30 

341  RCL  05 

366  X 

342  XZQ^BETA 

367  RCL  33 

343  X2QTG1 

368  / 

344  mfaz 

369  FIX  2 

345  RCL  32 

370  ^DEL  T 

346  RCL  28 

371  ARCL  X 

347  1/X 

372  AV1EW 

348  YtX 

373  RTS 

349  STO  33 

374  END 

350  RCL  28 

A. 


APPENDIX  C 


CALCULATION  OF  TEMPERATURE,  INTERNAL  ENERGY,  AND  ENTHALPY 
FOR  ISENTROPIC  FLOW 

In  ord«r  co  calculate  the  flow  field  downstream  of  the  shock  wave, 
as  illustrated  in  Figure  C-l,  one  needs  to  calculate  the  internal  energy 
and  enthalpy  of  the  water  as  it  expands.  The  energy  equation  appropriate 
for  the  flow  illustrated  in  Figure  C-l  is 

|j(pu(e  +  - — y^-)  +  pu]  +  |p[pv(e  +  - — y^-)  +  pv] 

2  2  (C-l) 

+  •ptpv(e  +  u  ^  )  +  pv]  ■  0 

Equation  (C-l) ,  which  assuaes  inviscid  flow  without  heat  conduction,  indicates 
the  need  to  calculate  Internal  energy. 


Figure  C-l.  Sketch  of  Flow  Field  Due  to  an  Aluminum  Jet  Penetrating  Water. 
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A.  Calculation  of  Internal  Energy  and  Enthalpy 

Sy  definition  internal  energy  and  enthalpy  are  related  by 

H  ■  e  +  pv  (C-2) 

Capital  H  la  uaad  for  enthalpy  since  h  la  used  for  dissipated  enthalpy. 

Refer  to  Figure  C-2  which  is  sin liar  to  Figure  1  of  Appendix  A.  The 
points  indicated  in  Figure  C-l  are  shown  also  in  Figure  C-2.  For 
example,  in  Figure  C-l  point  3  is  inmediataly  downstream  of  the  shock  wave. 

In  passing,  one  should  note  chat  to  calculate  shock  properties,  the 
normal  component  of  velocity  should  be  used.  Point  3  is  shown  in  Figure  C-2. 
As  water  flows  along  the  streaatube,  expansion  occurs.  The  pressure 
decreases  as  indicated  by  point  4  in  Figure  C-l  or  Figure  C-2. 


temperature. 


pressure,  p 


Figure  C-2.  Eugonlot  and  Mia be tic  (Isentropic)  Curves  in 
the  Pressure-Temperature  Plane. 
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From  equation  2.20  of  Appendix  A 


AH  •  H3 


H0  *  “31  +  h10 


(03) 


The  subscripts  31,  used  above  and  later  in  other  equations,  seen  jj 
is  evaluated  along  the  adiabatic  curve  between  points  3  and  1. 
Likewise  is  evaluated  along  a  line  of  constant  pressure  between 
points  1  and  0  in  Figure  02.  From  equation  (2.19)  of  Appendix  A 


u> 


31 


n-1 


-  1] 


(C-4) 


Also,  from  equation  (2.20)  of  Appendix  A, 

1 

h,Q  -  /  C  (0,T')dT* 
0  P 


(05) 


At  some  intermediate  point,  such  as  point  4  in  Figure  02,  one  replaces 
subscript  3  by  4  in  equations  (C-4)  and  (05). 

Combining  equations  (C-2)  through  (05),  an  expression  for  Internal 
energy  e  can  be  derived;  the  result  is 


2  a11"*1 

*  "  h10  +  Clln(n  -  1) 


n:  n-1 


(C-6) 


where  z  is  v^/v.  Hence  equation  (C-6)  uses  v  as  the  Independent  variable. 
An  alternate  form  can  be  derived  using  p  as  the  independent  variable;  the 
reeult  is 


n-1 

n 


•  -  ho  * 


,1/n  n-1 


(07) 


n(l  +  x)J 

where  x  ■  p/3  *  p/A.  The  symbols  3  and  A  are  defined  on  page  22,  Appendix  A 
see  equation  (XI-2).  Equations  (C-6)  and  (07)  have  a  very  similar  form. 


An  universal  equation  can  be  written  as 

2  w11"1  1  1 

*  "  **10  +  ClCn(n  -  1)  +  nw  "  n^T1 


(08) 


1 


where 


z  for  v  /v  as  independent  variable 
w  »  {  ..  (C-9) 

(1  +  x)  for  p/B  as  Independent  variable 

Note  thee  is  a  constant  value  for  any  streamtube. 

B.  Calculation  of  Temperature 

If  cha  energy  equation  incorporates  thermal  conduction,  a  term 
involving  temperature  is  used.  The  term  has  the  form 

Q  -  kfj  (C-10) 

2 

where  Q  is  the  heat  flux  in  Joule/m  sec,  and  k  is  thermal  conductivity. 

To  calculate  the  temperature  at  some  intermediate  point,  i.e.  T,, 
one  can  use  the  subroutine  DEL  T.  In  DEL  T,  p,  is  replaced  by  p^.  As 
a  note  of  caution,  when  intermediate  values  of  temperature  are  calculated, 
the  prasaura  behind  cha  shock  wave,  p^,  is  replaced  in  storage  register  12 
by  p4<  che  intermediate  pressure. 

C.  Subroutine  E-BAR  and  User's  Guide  to  E-BAR 

A  subroutine  has  been  written  to  calculate  Internal  energy  and 
enthalpy.  A  listing  of  the  subroutine  is  given  below: 


01 

lblte~bar 

11 

YtX 

21 

- 

02 

tisxn  p? 

12 

STO 

35 

22 

/ 

03 

PROMPT 

13 

RCL 

28 

23 

RCL 

28 

04 

STO 

36 

14 

1 

24 

1/X 

OS 

RCL 

08 

15 

- 

25 

RCL 

35 

06 

/ 

16 

Y+X 

26 

/ 

07 

1 

17 

RCL 

28 

27 

+ 

08 

+ 

18 

/ 

28 

RCL 

28 

09 

RCL 

28 

19 

RCL 

28 

29 

1 

10 

1/X 

20 

1 

30 

32 


31 

1/X 

41 

AVIEW 

51 

RCL 

35 

32 

- 

42 

ST0  37 

52 

/ 

33 

RCL 

17 

43 

STOP 

53 

RCL 

37 

34 

ENTER* 

44 

RCL  36 

54 

+ 

35 

X 

45 

1.013E8 

55 

th  bar 

36 

X 

46 

X 

56 

ARCL  X 

37 

RCL 

07 

47 

RCL  26 

57 

AVIEW 

38 

+ 

48 

0.001 

58 

ST0 

38 

39 

te-bar  - 

49 

X 

59 

END 

40 

ARCL  X 

50 

X 

As  a  convenience,  subroutine  E-BAR  should  be  assigned  to  some  key  when 
calculator  is  in  USER  mode. 

To  operate  E-BAR,  one  must  execute  HH  as  discussed  on  pages  4  to  9. 
Operation  of  HH  fills  the  registers  with  values  which  are  recalled  for 
use  by  E-BAR.  For  a  sample  problem,  HH  was  executed  using  n  *  7.15, 

P  *  90  kilo bars,  and  tQ  »  0.  When  HH  has  been  completed,  one  proceeds  as 
follows: 

1.  Continue  with  HP41CV  in  USER  mode. 

2.  Assign  E-BAR  to  a  key  on  the  keyboard  by  pressing 

1  1  ALPHA  E-BAR  ALPHA,  R+ 

In  this  example,  E-BAR  has  been  assigned  to  the  R+  key. 

3.  Press  R+,  and  observe  ISEN  ??.  This  is  the  isentropic  pressure 
at  point  4  in  Figure  C-2;  obviously  one  can  select  any  value. 

The  number  in  parentheses  is  the  value  used  in  the  sample 
problem  (60  kilobars) . 

4.  Press  R/S. 
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5.  Calculator  displays  E-BAR  **.  For  tha  sample  problem 
(E-BAR  -  1,249,766  Joule/kg) . 

6.  Press  R/S. 

7.  Calculator  displays  H  BAR  ■.  For  tha  example  problem 
(H  BAR  -  5, 449, 312  Joula/kg) . 

Subroutine  E-BAR  has  now  been  completed.  If  calculations  are  desired  at 
other  values  for  p^,  continue  as  follows: 

3.  Press  R/S,  and  observe  ISEN  ??. 

9.  Input  ISEN  P. 

10.  Press  R/S. 

11.  Etc. 

D.  Subroutine  ISEN  T  and  User's  Guide  to  ISEN  T 

A  subroutine  ISEN  T  has  been  written  to  calculate  temperature  along 
the  adiabatic  (is an tropic)  curve  of  Figure  C-2.  A  listing  of  the  subroutine 
is  given  below: 


01  lbltisen  t 

06  STO  39 

11  tisen  t 

02  TISEN  P? 

07  RCL  01 

12  ARCL  X 

03  PROMPT 

08  + 

13  AVTEW 

04  STO  12 

09  273.16 

14  END 

05  XEQTDEL  T 

10  + 

To  operate  ISEN  T,  one  must  execute  HH  as  discussed  on  pages  4-9.  For 
a  sample  problem,  HH  has  been  executed  with  n  •  7.15,  p  »  90  kilobars,  and 
tQ  -  0. 

1.  Continue  with  HP41CV  in  USER  mode. 

2.  Assign  ISSN  T  to  a  key  on  the  keyboard  by  pressing 

(“|  ALPHA  ISEN  T  ALPHA,  EEX 
la  this  example,  ISEN  T  has  been  assigned  to  EEX  key. 


3.  Press  EEX,  and  observe  ISEN  ?Y.  This  is  che  isentropic 
pressure  ac  point  4  in  Figure  C-2.  Since  p^  or  ISEN  P  is  an 
independent  variable,  any  value  can  be  selected.  The  number  in 
parentheses  is  Che  number  used  in  che  sample  problem  (60  kilobars) . 

4.  Press  E./S. 

5.  Calculator  displays  ISEN  T  ■.  For  the  sample  problem,  the  result 
is  (ISEN  T  »  740.6).  The  dimensions  are  °K. 

Subroutine  ISSN  T  has  nov  been  completed.  If  calculations  are  desired 
at  other  values  for  p^,  continue  as  follows: 

6.  Press  R/S,  and  observe  ISEN  P?. 

7.  Input  ISEN  P  (50  kilobars). 

3.  Press  R/S. 

9.  Observe  ISEN  T  -  691. 7°K. 

10.  Etc. 

E.  Additional  Assignment  of  Storage  Registers  for  Subroutines  E-BAR  and  ISEN  T 

Table  C-I  gives  the  additional  assignment  of  storage  registers. 


Table  C-I.  Additional  Storage  Registers 

for  E-BAR  and  ISEN  T 

Register 

Symbol 

Definition 

Units 

Programs 

35 

(1  +  p/B)1/a 

Also  symbol  w  in 
equation  (C-8) . 
Equals  v^/v. 

- 

E-BAR 

36 

P 

Intermediate 
pressure  on 
Isentropic  line. 

kilobars 

E-BAR 

37 

e 

Internal  energy. 

Joule/kg 

E-BAR 

38 

H 

Enthalpy. 

Joule/kg 

E-BAR 

39  AT  (Isentropic) 

Change  in  tempera¬ 
ture  along  the 
adiabatic  line  in 
Figure  C-2. 

ISEN  T 

35 


F.  Sample  Results  for  Properties  Along  Adiabatic  (Isen tropic)  Line 

Table  C-II  suamarizes  the  results  for  conditions  along  the  adiabatic 
line  in  Figure  C-II.  Input  values  for  the  shock  wave  were  n  ■  7.15,  p  *  90 
kilobars,  tg  »  0°C,  and  t^  ■  167.762°C. 


Table  C-II.  Sample  Results  for  Conditions  Along 
Adiabatic  (Isen tropic)  Line 


p 

kilobars 

e 

kilojoule/kg 

a 

kilojoule/kg 

o 

5*  H 

0 

674 

674 

441 

5 

707 

1180 

460 

10 

755 

1634 

486 

20 

857 

2474 

539 

30 

958 

3260 

591 

40 

1057 

4013 

642 

50 

1154 

4741 

692 

60 

1249 

5449 

741 

70 

1343 

6142 

789 

80 

1435 

6820 

836 

90 

1526 

7488 

883 

For  p  *  0  kilobars,  the  temp  jrature  T  is  the  residual  due  to  entropy  increase 
across  the  shock  wave.  For  p  ■  90  kilobars,  the  temperature  T  is  the  value 

i 

lamed  lately  behind  the  shock  wave. 
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As  the  needs  for  computational  flexibility  have  changed,  the  report  of 
Reference  (a)  has  grown.  Enclosed  is  Appendix  D,  which  permits  calculation 
of  pressure  as  a  function  of  specific  volume  and  Internal  energy.  An 
iteration  is  required,  and  an  efficient  Iteration  procedure  has  been 
developed.  A  fractional  error  for  the  iteration  variable  of  5E-05  can 
be  achieved  in  2  or  3  iterations. 
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APPENDIX  D 
(NPS67-82-001) 

CALCULATION  OF  PRESSURE  AS  A  FUNCTION  OF  SPECIFIC  VOLUME  AND  INTERNAL  ENERGY 
Introduction 

In  order  to  use  efficiently  the  computer  program  developed  by  Kutler  and 
Chakravarthy  [D-l],  the  ability  to  calculate  pressure,  p,  as  a  function  of 
specific  volume,  v,  and  internal  energy,  e,  is  needed.  The  program  for 
HP41CV,  which  is  described  in  this  appendix,  calculates  p  ■  p(v,e). 

Relevant  Equations 

The  relevant  equations  for  the  calculation  can  be  assembled  from  the  equations 
in  this  report.  Define 

z  ■  Vj^/v  (D-l) 

where  v^  is  the  specific  volume  at  t^,  and  v  is  the  specific  volume  at  the 
desired  condition.  The  value  for  v^  can  be  calculated  from 

v(t)  -  0.99415  +  0. 0002929 (t  -  25)  +  3.241  x  10-6(t  -  25)2  (D-2) 

The  internal  energy  is  given  by 


—  .  .2.  z 

e  -  hio  +  c, l 


n~l  11  — 

,  i  -/  : - rr  + - r]  *  e(t.  ,v.  /v) 

1  n(n  -  1)  nz  n  -  1  11 


(D-3) 


Note  that  e  is  a  function  of  t^  and  v^/v.  The  quantity  h^Q  is 


h10  "  f  C  <0.t)dt  -  3.9644(tx  -  tQ)  +  3.12  x  10_4(t*  -  t2) 
fc0 


(0-4) 


The  temperature,  tg,  is  the  reference  temperature  which  can  be  the  initial 
water  temperature  before  the  shock  wave  arrives. 

The  pressure  can  be  related  to  specific  volume  by  the  Tait  equation 

p  -  BUjHCVj/v)"  -  1]  -  p(v,tx)  (D-5) 
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r 


The  parameter  B(t)  is  given  by 

B(t)  -  3.1' 4  -  1.65  x  10~3(t  -  55)  -  1.181  x  10~4(t  -  55)2 

7  ,  (D-6) 

+  5.32  x  10  '(t  -  55) 

2 

Finally  for  equation  (D-3),  the  value  of  c^  must  be  calculated.  The 
appropriate  formula  is 

C1  "  nviB^i^  (D-7) 

The  value  of  n  is  7.15. 

Iteration  Procedure 

The  preceding  equations  permit  calculation  of  p  =  p(e,v); 

however,  an  explicit  formulation  is  not  possible,  or,  at  least,  is  very 

difficult.  Consequently,  an  Iterative  approach  is  used.  Define  e^  as  a 

calculated  value  from  equation  (D-3).  Define  e^  as  an  input  value  at  which 

the  pressure  is  desired.  As  implied  by  Figure  D-l,  an  iterative  procedure  is 

needed  to  find  the  correct  value  for  t,  at  which  e  *  e , . 

c  i 

The  function  e \  is  a  monotonically  increasing  function  in  the  temperature 

region  of  interest.  An  Iteration  scheme  using  straight  lines  is  used  to  find 
t^.  For  computer  programming,  capital  letters  are  used;  hence,  in  terms  of 
capital  letters,  the  equation  for  the  straight  line  is 

T3,J  -  f^rfjCT2  -  Tl)  +  T1  (D-8) 

where  the  symbol  definitions  are  as  follows: 

T3,J  calculated  value  for  temperature  using  a  straight 
line  as  shown  in  Figure  D-2;  J  is  an  index  which 
indicates  the  number  of  the  iteration. 

E  desired  value  of  E  BAR;  also  identified  as  input  value  e^. 

El  value  of  E  BAR  at  temperature  Tl. 

E2  value  of  E  BAR  at  temperature  T2. 
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Figure  D-l.  Iteration  to  determine  correct  value  of  t 


(a)  Initial  straight  line 


(b)  Calculation  of  E3,l 


(c)  Revised  straight  line  (d)  Calculation  of  E3,2 

and  solution  for  T3,2 


Figure  D-2.  Illustration  of  iteration  procedure. 


I 


Refer  to  Figure  D-2(a).  The  first  values  for  tl  and  T2  are  inputs  and  are 
guesses.  The  program  calculates  El  and  E2  and  solves  for  T3>1  as  shown  in 
Figure  D-2(a).  The  value  for  T3,l  is  used  to  calculate  E3.  E3  is  compared 
to  E.  If  the  difference,  E  -  E3,  is  too  large,  then  T2  «  T3,l,  Tl  *  T2, 

E2  =  E3,  and  El  ■  E2.  A  new  straight  line  is  calculated  passing  through 
points  b  and  d  in  Figure  D-2(b)  and  Figure  D-2(c).  A  solution  is  obtained 
for  T3,2,  which  is  the  second  value  for  T3.  Using  T3,2,  E3  is  calculated  and 
is  once  again  compared  with  E.  If  the  quantity  ABS(E  -  E3)  is  too  large,  the 
values  for  E  and  T  are  reassigned  as  follows:  T2  -  T3,2,  Tl  *  T2,  E2  =  E3, 
and  El  -  E2. 

Once  again  a  solution  is  obtained  for  T3,3  using  a  line  passing  through 
points  d  and  e  in  Figure  D-2(d) .  The  procedure  repeats  until  the  desired  accuracy 
for  ABS(E  -  E3)  is  obtained. 

Program  Description 

Figure  D-3  is  a  flow  chart  for  the  program  p  -  p(e,v),  which  has  been 
labelled  P  EV.  Additional  storage  registers  are  assigned  as  shown  in 
Table  D-l.  The  assignment  of  registers  1  to  39  is  given  in  Table  1  and 
Table  C-l.  Table  D-2  is  a  program  listing.  The  various  subroutines,  such  as 
VW,  BT,  and  H10,  must  be  in  the  memory  of  HP41CV. 

As  shown  in  the  flow  chart  of  Figure  D-3,  inputs  are  made  for  B,  E  BAR, 

J,  Tl,  and  T2.  An  input  J  ■  0  must  always  be  used;  otherwise  the  program  logic 
will  be  in  disarray.  The  quantities  T  and  V  are  the  variables  used  in  subroutine 
LBL  12  to  calculate  E.  The  quantities  utilized  or  calculated  as  a  function 


of  J  are  as  follows: 

VALUE  OF  J 

VALUE  OF  T 

VALUE  OF  E 

0 

Tl 

El 

1 

T2 

E2 

2 

T3 

E3 
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Table  D-l.  Additional  Storage  Registers 

for  P  EV 

Register 

Symbol  Definition 

Units 

Programs 

40 

V 

Specific  volume  at  which  pressure  is 
desired 

cm  /  gm 

12 

41 

E 

Internal  energy  at  which  pressure  is 
desired 

kilojoule/kg 

12 

42 

T 

Value  of  t,  as  used  in  equation  (C-5) 

°C 

VW, 

BT 

43 

- 

Not  assigned 

- 

- 

44 

V1 

Value  of  specific  volume  at  p  ■  0  and 
temperature  t^ 

cm'Vgm 

12, 

13 

45 

B(t) 

Parameter  in  equation  of  state 
relating  p  to  v 

kilobars 

12, 

13 

46 

2 

C1 

Speed  of  sound  squared  at  tempera ture 
*1 

2/  2 
m  /sec 

12 

47 

h10 

Contribution  to  enthalpy  at  zero 
pressure  between  temperatures 

to  and  h 

kilojoule/kg 

12 

48 

z 

Ratio  v^/v 

- 

12 

49 

El 

Internal  energy  evaluated  at  T1 

Joule/kg 

12 

50 

E2 

Internal  energy  evaluated  at  T2 

Joule/kg 

12 

51 

T1 

Temperature  at  point  1  on  straight 
line 

°C 

12 

52 

T2 

Temperature  at  point  2  on  straight 
line 

°c 

12 

53 

J 

Index  to  count  number  of  executions 
of  LBL  12 

- 

09 

54 

E 

Register  to  store  output  of  LBL  12 
which  is  E  ■  E(T,v) 

Joule/kg 

12 

55 

T3 

Temperature  of  point  determined  by 
straight  line 

°C 

11 

56 

E3 

Internal  energy  at  temperature  T3 

Joule/kg 

12, 

14 

46 


Table  D-2.  Program  Listing 


01 

lbltp  ev 

26 

xeqtbt 

51 

1/X 

76 

RCL 

53 

101 

RCL 

51 

126 

RCL  50 

02 

T 

V? 

27 

STO  45 

52 

7.15 

77 

1 

102 

- 

127 

STO  49 

03 

PROMPT 

28 

7.15 

53 

/ 

78 

+ 

103 

X 

128 

RCL  56 

04 

STO  40 

29 

X 

54 

+ 

79 

STO 

53 

104 

RCL 

51 

129 

STO  50 

05 

te  BAR? 

30 

RCL  44 

55 

6.15 

80 

GTO 

11 

105 

+ 

130 

GTO  11 

06 

PROMPT 

31 

X 

56 

1/X 

81 

LBL 

09 

106 

T 

T3 

= 

131 

LBL  13 

07 

1000 

32 

1000 

57 

- 

82 

RCL 

54 

107 

ARCl 

,  X 

132 

7.15 

08 

X 

33 

/ 

58 

RCL  46 

83 

STO 

49 

108 

AVIEW 

133 

RCL  48 

09 

STO  41 

34 

STO  46 

59 

X 

84 

RCL 

52 

109 

STOP 

134 

XOY 

10 

T 

J? 

35 

RCL  42 

60 

RCL  47 

85 

STO 

42 

110 

STO 

55 

135 

Y-t-X 

11 

PROMPT 

36 

xeqth  10 

61 

+ 

86 

RCL 

53 

111 

STO 

42 

136 

1 

12 

STO  53 

37 

1000 

62 

te  = 

87 

1 

112 

GTO 

12 

137 

- 

13 

T 

Tl? 

38 

X 

63 

ARCL  X 

88 

+ 

113 

LBL 

14 

138 

RCL  45 

14 

PROMPT 

39 

STO  47 

64 

AVIEW 

89 

STO 

53 

114 

RCL 

54 

139 

X 

15 

STO  51 

40 

RCL  44 

65 

STOP 

90 

CTO 

12 

115 

STO 

56 

140 

1.013E8 

16 

T 

T2? 

41 

RCL  40 

66 

STO  54 

91 

LBL 

11 

116 

RCL 

41 

141 

/ 

17 

PROMPT 

42 

/ 

67 

RCL  53 

92 

RCL 

50 

117 

- 

142 

TPRESS  ' 

18 

STO  52 

43 

STO  48 

68 

X-0? 

93 

RCL 

49 

118 

ABS 

143 

ARCL  X 

19 

RCL  51 

44 

6.15 

69 

GTO  09 

94 

- 

119 

50 

144 

AVIEW 

20 

STO  42 

45 

YtX 

70 

1 

95 

RCL 

41 

120 

X>Y? 

145 

END 

21 

LBL  12 

46 

7.15 

71 

RCL  53 
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Once  J  has  acquired  a  value  of  2,  the  program  loops  repeatedly  as  follows: 

LBL  12,  TEST  ABS(E-E3),  LBL  14,  LBL  11,  LBL  12,  etc.  When  the  criterion 
for  the  value  of  ABS(E-E3)  is  attained,  the  program  departs  from  the  loop 
and  moves  to  LBL  13. 

Sample  Run 

The  program  P  EV  has  been  assigned  a  key  in  USER  mode  of  HP41CV.  Pushing 
that  key  causes  the  following  queries  to  appear: 


V? 

0.75  cm^/gm 

E  BAR? 

900  kilojoule/kg 

J? 

0 

Tl? 

50°C 

T2? 

150°C 

Once  the  input  has  been  made,  the  crow’s  foot  moves  across  the  display  indicating 
the  program  is  being  executed.  Next  to  appear,  after  R/S  is  pushed,  is 

E  *  377,507  (Joule/kg) 

which  is  the  value  for  El  calculated  at  Tl.  Next 

E  ■  880,354  (Joule/kg) 

appears  which  is  E2.  The  situation  for  the  straight  line  is  exactly  as  shown 
in  Figure  D-2(a) .  Pushing  R/S  yields 

T3  -  153.91  («  T3,1°C) 

P-’hing  R/S  causes  the  program  to  display 

E  -  902,126 

which  is  E3.  The  situation  is  that  illustrated  in  Figure  D-2(b).  The  sequence 
of  values  is  listed  below 

T3  -  153.53  (-  T3,2°C) 

E  -  899,985 

PRESS  -  30.45  (kilobars) 
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Reference  to  line  119  of  Program  Listing  indicates  that 

jE-E3|<50 

has  been  achieved. 

To  indicate  the  influence  of  the  initial  values  for  Tl  and  T2,  the 
program  was  run  again  with 

V  «  0.75  cm^/gm 
E  BAR  =  900  kilojoule/kg 
J  =  0 
Tl  =*  0 
T2  =  10 

The  sequence  of  values  is  listed  below: 

E  *  135,345  (=  El) 

E  -  184,447  («  E2) 

T3  *  155.73  (=  T3,l) 

E  -  912,390  («  E3) 

T3  <=  153.25  (**  T3.2) 

E  =  898,425 
T3  •  153.53  («  T3,3) 

E  -  899,993 
PRESS  =  30.45 

Using  values  of  Tl  and  T2  quite  far  removed  from  final  value  of  T3,J,  the 
program  converged  to  a  solution  in  three  iterations. 
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